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I. INTRODUCTION 



Quite recently there has been a notable in- 
crease of experimental activities aiming at the 
exploration of properties of thermally gener- 
ated fluctuating electromagnetic fields close to 
the surface of some material, and at detect- 
ing the near-field mediated heat transfer ii~— 
Hu et al. have measured the near-field ther- 
mal radiation between um-spaced glass plates, 
and have demonstrated that the resulting near- 
field heat transfer exceeds the far-field limit set 
by Planck's blackbody radiation law.— Next, 
Narayanaswamy et al. and Shen et al. have 
studied the heat transfer between microspheres 
and flat substrates, with emphasis on the cou- 
pling of surface phonon polaritons across the 
gap between them, and have reported heat 
transfer coefficients three orders of magnitude 
above the blackbody radiation limit Then 
Rousseau et al. have carried out precise mea- 
surements of the radiative heat transfer be- 
tween sodalime glass spheres with diameters of 
22 or 40 um and borosilicatc glass plates for 
distances ranging from 30 nm to 2.5 um,— and 
have verified theoretical predictions based on 
fluctuational electrodynamics^ with impres- 
sive accuracy. On the other hand, significant 
progress has been made at using near-field ef- 
fects for imaging. Kittel et al. are develop- 
ing a device termed Near-Field Scanning Ther- 
mal Microscope (NSThM)2r— which does not 
yet seem capable of highly accurate quantita- 
tive measurements of the near-field heat cur- 
rent between its sensor and the sample, but 
which lends itself to nanoscale thermal imaging 
of structured surfaces.— Moreover, De Wilde 
et al. have reported the successful operation 
of a Thermal Radiation Scanning Tunneling 
Microscope (TRSTM),^ providing images of 
thermally excited surface plasmons, and giving 
clear evidence for spatial coherence effects in 



near-field thermal emission. 

These remarkable developments indicate that 
thermal near-field physics, after having been 
under intense theoretical investigation for some 
time already,— ~—£ is breaking through to the 
forefront of experimental research right now. 
There are several compelling reasons for this 
trend: Besides the prospects of obtaining novel 
insight into fundamental physics in dielectric 
matter, and of developing advanced diagnostic 
tools for materials science, thermal near-field ef- 
fects have great potential for near-field thermo- 
photovoltaic energy conversion^ - — 

On the theoretical side, one of the most im- 
portant quantities characterizing the fluctuat- 
ing thermal near field close to a dielectric sur- 
face is its local density of states (LDOS)j^ In 
particular, the power P transferred between 
a dielectric sample at temperature Tg and a 
nanoparticle at temperature Tp, brought into 
the sample's near field at a position a such that 
the particle may effectively be treated within 
the dipolc approximation, and the heat transfer 
proceeds almost entirely via evanescent modes, 
is given by (see, e.g., Refs. [20l - [25h 



P = / duj 2u [Q(u, T P ) - 0(w, T s )] (1) 



x[4( W )JJ E ( W) a)|^Hfl ! ( W ,a)] , 

where D E (w,a) is the electric and Z) H (w,a) 
is the magnetic part of the sample's LDOS at 
the point a of effective interaction; «p(w) and 
/ip(w) denote the imaginary part of the par- 
ticle's electric and magnetic polarizability, re- 
spectively. Finally, 



e(w,r) 



exp(huj I ksT) — 1 



(2) 



is the Bose-Einstcin function; the sign in Eq. (QJ 
is chosen such that a net energy transfer from 
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the particle to the sample, occurring for Tp > 
Ts, gives a positive P. The use of the dipole 
approximation underlying Eq. (QJ requires that 
the distance of the nanoparticle from the sam- 
ple's surface remains large compared to its lin- 
ear size. Those frequencies which significantly 
contribute to the heat transfer are limited by 
the higher temperature T max = max(Tp, Tg). 
Provided the polarizabilitics and the LDOS ex- 
hibit no resonances in the accessible frequency 
regime, requiring in particular the absence of 
thermally excitable surface modes, the main 
contribution to the integral ([T]) merely stems 
from frequencies in the vicinity of the thermal 
frequency Wth ~ 2.82 fceT/fi, so that 

P gc uj th [a'p(uJth)D E {uj th ,a) (3) 
+ ^(cj th )D H (w th ,a)] . 

Under suitable conditions, already this simple 
approximation can give surprisingly good re- 
sults when trying to theoretically reconstruct 
surface images obtained with the NSThMJ^ 

The possibility to experimentally assess the 
LDOS above nanostructured surfaces demands 
refined techniques for its calculation. Assum- 
ing local thermal equilibrium, and consider- 
ing positions r so close to the sample's surface 
that the energy density is dominated by evanes- 
cent modes and the contribution of propagating 
modes can be neglected, the electric and the 
magnetic LDOS are related to the imaginary 
parts of the traces of the renormalized (or 're- 
flected') electric and magnetic Green's dyadics 
Gf and Gj 3 through the relations^ 

D e (lj, r) = -^-ImTrGf (r,r) (4) 

7TC Z 

and 

D K (co,r) = -^■ImTrG* I (r,r) . (5) 

7TC 

In the present paper we exploit this connec- 
tion for computing the LDOS above a nano- 
structured surface to second order in the surface 
profile, relying on an earlier formulation of the 
perturbation scries by Grcffct.— The method 
is technically involved, and soon hits practical 
computational limits when proceeding to higher 
orders. Nonetheless, we show that second-order 
calculations now are feasible routinely. Our 
work thus extends a previous study,— which has 
given first-order results, and enables us to delin- 
eate under which conditions low-order pertur- 
bation theory is sufficient. It also complements 
a recent investigation by Biehs and Grcffct 28 
who have considered rough surfaces described 
as stochastic Gaussian processes. In contrast, 



we focus on surfaces with deliberately induced, 
regular nanoprofilcs. We proceed as follows: In 
Sec. |n] we sketch the underlying perturbativc 
scheme^ and outline a few details of its nu- 
merical implementation, deferring technicalities 
to the Appendix [A] We then present results 
of our computations in Sec. IIII1 first examin- 
ing the relative magnitude of first-and second- 
order contributions, and then outlining how to 
quantify the resolution power of an idealized 
NSThM. Some conclusions are drawn in the fi- 
nal Sec. HVl 

II. COMPUTATION OF THE GREEN'S 
DYADICS 

In this section we utilize an analytical pertur- 
bative approach, originally developed by Gref- 
fet for calculating the scattered electromagnetic 
waves above a rough dielectric surface,— in or- 
der to obtain the required Green's dyadics Gf 
and Gj 3 . Grcffct's approach results in a series 
of recursively determined contributions in as- 
cending orders of the surface profile, and thus 
enables one to systematically assess the higher- 
order terms. 



A. Calculational scheme 

We assume that the surface is described by an 
expression z ~ hf(x, y), where /(x, y) is a nor- 
malized function varying between +1 and — 1; 
the scaling parameter h carries the dimension of 
a length. The nonmagnetic dielectric medium, 
equipped with permittivity e(w), fills the entire 
half-space z < hf(x,y). For z > hf(x,y), out- 
side the dielectric, the total electric field con- 
sists of a prescribed incident component E\(r). 
and of the so far unknown reflected compo- 
nent E r (r), while the transmitted field inside 
the dielectric is denoted as E t (r). Greffet has 
given a recursive series solution for the trans- 
mitted and the reflected field^ invoking the 
extinction theorem and the Raylcigh hypoth- 
esis.— The extinction theorem amounts to an 
exact integral formulation of the boundary con- 
dition, whereas the use of the Rayleigh hypoth- 
esis means expanding the transmitted, incident, 
and reflected fields in plane waves travelling in 
z-direction, 

E t (r) = U 2 Ke t {K)e^P- k ^ , (6) 



Ei(r) = / d 2 Kei (K)c i{K - p - k * oz) , (7) 
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and 



/ 



E T {r) = / d 2 /te r (K)e 



(8) 



Here we write r = [x,y, z] for the position 
vector, p = [x, y, 0]* for its lateral part, and 
K = [kxjkyjO]* for the lateral component of the 
wave vector; moreover, 



k z = y/ eko - K 



and 



(9) 



(10) 



ei(/c) = e i)S (fv)a s (rc) + e iiP (K)a pi0 («) , (15) 
and 

e r (K) = e r , s (K)a s (rc) + e riP (rc)a+ (rc) , (16) 



where 



a s («) 



o p ,t( K ) 



Mi- 




kyk z 



(17) 



(18) 



are the ^-components of the wave vector in- 
side and outside the medium. We also use the 
notation fco = w / c an d K = 1^1- These ex- 
pansions ([HI)-© assume translational symme- 
try in the x-y-plane and hence are strictly jus- 
tified outside the structured region, that is, for 
z > h and z < —h. However, the extinction 
theorem requires to evaluate the fields on the 
very surface of the dielectric, where the valid- 
ity of the above expressions cannot be taken 
for granted. Ignoring this complication and 
using the expansions ©-(HI) nonetheless is a 
common procedur e) 30 ' 31 which has been looked 
into by several authors from the mathemati- 
cal point of viewj 3 ^— it appears to work re- 
liably at least for sufficiently small values of h. 
For example, in the case of a sinusoidal grat- 
ing described by z — (h/2) cos(2irx/D) this 
Rayleigh hypothesis holds rigorously 3 ^— up to 
h max /D = 0.142 521, and may therefore be em- 
ployed for both the propagating and the evanes- 
cent parts of the field as long as the ratio h/D 
stays below this boundary. 

The field's Fourier components then are ex- 
panded in the forms 



e t («) 



m=0 



Am) 



(«) 



(ii) 



(12) 



and 



e i( K ) = E 



m=0 



(m) i ■ 

e> '(h) _ ( o) 



rw- (13) 



It is useful to split the fields into their s- and 
p - components according to 

e t (rc) = e ti8 (/c)a s (K) + e t , p (K)a~ t (K) , (14) 



a P ,o(«0 = - 



1 



and 



1 



k x k z o 
kyk z o 



k x k z o 
kyk Z Q 



(19) 



(20) 



here n = y/e is the index of refraction. 

Following Grcffct,^ one then obtains the 
transmitted field in the recursive form 



e t ( ?(") 



kz °R- 1 (K,K) / dV <J R(«, «') 



4tt 2 

'm\ - k) 



(m-q) i , 



-t.s 



(*') 

a (rn-q) f^i 
-t,p 



(21) 



so that e[ (k) is proportional to h m ; the ac- 
cording expression for the reflected field reads 



e^ } («) 



8^{/ dV (^'^ 



*E 

eS(K') 
+ P(k,k)- 



m\ (-i/i) m -9F( m -«)(K' - re) 
q) (k zQ + KY +q ~ m 



4tt* 



'20 



Am) 
-t.s 
Am) 



(«) 
p» 



(22) 



In these equations the quantities F( n \n) de- 
note the Fourier transforms of powers of the 
surface function, 



fW((t) = J d 2 P c 1K - p f n {p) 



(23) 
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The linear operators H(k, k') and P(k, k') 
effectuate the double vectorial product with 
fe r = [k x ,ky, fc z o]' and fc r = [k Xl ky, fc z o]*, 
respectively; these double products (namely, 
k~ x k~ x and k^xk^x) typically appear when 
using the extinction theorem. In the basis cho- 
sen, the matrix forms of these operators are 



R(k,k') 



and 



P(k,k') = 



Z-2 



/« • K 

fe^o (re x k/ 
fen 



fe, (re x re ) s 
nfeo 

k 2 k' 2 +n- n' k z ok' z 



nfe- 



(24) 



2/2 / ;,, ? / 

K K — f£ ■ rv - (") a' _ 

P " 



i ik 



(25) 

writing (k x k') z for the z-component of the 
vectorial product of k and k' . 

The above recursions start from the well 
known half-space results obtained for a per- 
fectly flat surface, which can be cast into the 
forms 



't.s 

-t,p 



and 



(«) rt s (K) o 

(k) L 



r,(«) 
r p («) 



ei, 8 (K) 

e i,p( K ) 



(26) 



ei, s (K) 



with the Fresnel coefficients 
2fc 20 



* s (k) 
r s (/c) = 



k z0 



k z 



From these fields 



*p(k) 



>(«) = 



2nfc g0 
n 2 fc z0 - 
n 2 fc z0 



(27) 



n 2 k z0 



?T|) -([22 |) we now proceed 
to the calculation of the Green's dyadics. More 
precisely, in order to compute the local density 
of states we need to determine the 'reflected' 
part of the Green's dyadics for coinciding source 
and observation points. To this end, we take 
the field of a delta-like source current located 
at p + ze z and pointing into the direction of 
the unit vector j as incident field, giving 



ei, s («) 

e i,p( K ) 



= c i(n-p+k z0 z) 



2k 



zO 



a p,o( K ) ' 3 



(28) 

To zeroth order in h, the reflected field then is 

(29) 



e$(«) 



^0_ c i(n-p+k z0 z) 

2kzo 



r B (K) 
r p («) 



a p,o( K ) ' 3 



The directions of the s- and p -components 
of the incident field are given by a s (re) and 
q~ (k), while the components of the reflected 
field are given by a s (/c) and q (k). Therefore 
it is convenient to split the Fourier coefficients 
of the Green's dyadics into the four parts that 
result from taking the dyadic products of these 
unit vectors, leading to 

9[ E '(»)( K )= Sl E i»'( K )a s W®a s W (30) 
+ g? : J£\K)a s (K)®a p0 (K) 

+ .9r E p ( s 0) («)ap,o(«)® a S («) 
+ g?j°\ K )a+ (K) <g> o- (k) . 

Having employed a delta-like source current 
the relation between the field and the electric 
Green's dyadic simply reads G E • j = E/ (iw/zo), 
so that the coefficients of this dyadic can eas- 
ily be read off from Eq. (|2"9"|) . By means of an 
inverse Fourier transform, equating source and 
observation point, one then arrives at the famil- 
iar result for the reflected Green's dyadic per- 
taining to a flat surface, 



G^> = -^/d 2 . e 



2ife~ 



2k z0 

r p( K )^o( K )® a p,o(' t )] 



[r s (/c)a s (fc) <x> a s (/c) 
(31) 



With the help of Eq. (|22j) one obtains similar 
expressions to all orders in the profile height h. 
An advantage of this approach lies in the fact 
that it is then quite easy to also determine the 
magnetic Green's dyadic, which is related to the 
electric one through^ 



G?(r,r) 



f 



2 VxGf(r,r')xVV=r • (32) 



In Fourier space the operator V is replaced ei- 
ther by i(« + k z0 e z ) or by i(« - k z0 e z ), de- 
pending on whether the curl acts on a unit vec- 
tor belonging to the incident or to the reflected 
field. Therefore, the magnetic Green's dyadic 
is derived from its electric counterpart by sim- 
ply replacing a(n) <S> a(n') by the expression 
pr(K + k z0 e z ) x o(k) ® o(k') x (k' - k' z0 e z ). 

Here we introduce k' and k' z0 because it is only 
to zeroth order that the reflected and the inci- 
dent field are characterized by the same wave 
vector. 

Finally, for calculating the trace of the 
Green's dyadics one just has to replace the 
dyadic products by their traces. Hence, the 
method sketched here yields a transparent 
strategy for obtaining the electric and the mag- 
netic LDOS above a structured surface. In 
practice, a restriction on the maximum order 
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achievable is imposed by the available compu- 
tational resources: To zeroth order only two- 
dimensional integrals have to be evaluated, to 
first order four-dimensional integrals appear; to 
second order one already has to deal with six- 
dimensional integrals, and so on. Clearly, a 
good choice of the numerical tools is decisive 
here; we therefore present some details of our 
implementation. 



B. Numerical implementation 

We exemplarily discuss the calculation of the 
electric local density of states to first and sec- 
ond order; its magnetic counterpart, and the 
higher-order terms, are determined in a similar 
way. With the help of Eq. (f22j) the trace of the 
first-order contribution to the electric reflected 
Green's dyadic is computed as 

Tr{G^ 1 )} = -|0e i ^FW(- g )a 1 (g) J 

(33) 

where a±(q) is given by the integral 

Oi (q) = J d V (q + K ; , k') ■ (g + «') 

(34) 

with the four-dimensional vectors S^ 1 ' and 
A^' specified in Appendix [X] We either em- 
ploy an experimentally determined surface pro- 
file,— or some suitably selected model function 
/(p); sample it, and perform a discrete FFT in 
order to determine . Then we compute a(q) 
for each required q by numerical integration, 
and finally take an inverse FFT of F^> (-q)a(q) 
to get the trace (f3"3"l) of the first-order Green's 
dyadic. 

To second order one has to deal with two con- 
tributions, one containing F@\ the other feed- 
ing from two factors F^ 1 ' . The former contri- 
bution has the same structure as the first-order 
term, 

Tr{G^ 2) }=|0e^F( 2 )(- g )a 2(g ), (35) 



with 



a 2 {q) = /dVS^q+K',K'Mr(gW,K') , 

(36) 

and is evaluated in the same manner. The other 
contribution contains a further integral, 

^ 2) }= f^^H-q) \% 



(E), 



(37) 



with integration variables q = k — k' and q' = 
k' — k" , and with the expression 

as (q, q') = J dV' S k ^+K a )- ( 38 ) 

X S 2 (^i Ki j Ki ) ■ -^4.^.^. (Av, ) , 
(2) (2) 

again, S\ and S 2 are stated explicitly in Ap- 
pendix El After computing a 3 (q,q') on a dis- 
crete mesh of q and q', a four-dimensional in- 
verse FFT is performed for determining the re- 
sulting contribution to the trace of the Green's 
dyadic. 

The integrals are numerically evaluated using 
Cuba routines^ the Fourier transforms are ex- 
ecuted by means of the FFTW library.— 

III. RESULTS 

In this section we discuss some numerical re- 
sults for the local density of states above exam- 
ple topographies, calculated to second order in 
the profile height. We first consider the rela- 
tive magnitude of the individual contributions, 
in order to estimate under which conditions the 
termination of the perturbation series can be 
justified. We then use the second-order data 
to discuss the resolution power of an idealized 
NSThM in a mode of operation in which the to- 
tal heat transfer is kept constant while scanning 
a sample's surface. 



Magnitude of second-order 
contributions 



Our basic model structure is a bar with 
height h and smoothed edges placed on an 
otherwise perfectly planar surface, infinitely 
extended in y-direction and possessing the 
width w in ^-direction, as described by the 
function 



hfi(x) 



1 



exp [C(|a 



0.5iy) 



1 



(39) 



Our calculations are done for h = 5 nm, 
w = 30 nm, and inverse smoothing length £ = 
10 9 m _1 . For comparison, we also consider the 
somewhat rounder profile 



hfi (x) = h exp 



1 



1 - {x/vf 



1 



(40) 



c ^PF^(-q')a 3 (q,q'), 



with v adjusted such that the respective areas 
under the two functions (f3"9"|) and (|40[) coincide. 
The resulting profile shapes are drawn in Fig. [TJ 
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Figure 1: (Color online) Profile function (|39p with 
parameters as employed in our model calculations, 
height h = 5 nm, width w — 30 nm, and inverse 
smoothing length £ = 10 9 m _1 (full line), together 
with the reference profile (|40[) (dashed). In either 
case, the dielectric properties of the sample are 
given by the Drude model with parameters for gold 
at 300 K. 



Because these profiles depend on only one 
variable, the Fourier transforms F^ n > contain 
delta functions, so that the integrals over q and 
q' in Eqs. ([55]). (f55j) . and (JHTJ) become effec- 
tively one-dimensional, drastically reducing the 
numerical effort. The profiles are discretized 
with a stepsize of 1 nm, covering a total range 
of 500 nm; we have checked that the numerical 
results thus obtained are stable against further 
reduction of the grid size to 0.5 nm. We assume 
that the dielectric function e(ui) of the samples 
is given by the Drude model^ 



(41) 



with plasma frequency w p = 1.4 x 10 s and 
inverse relaxation time 7 = 3.3 x 10 13 s _1 , de- 
scribing gold at the temperature T = 300 K. 

In Fig. [5] we depict the zcroth-, first-, and 
second-order contributions to the electric and 
to the magnetic part of the LDOS for the struc- 
ture ([STj]) . evaluated at a constant height of 
10 nm above the base plane for cj = 10 14 s _1 , 
roughly equal to the dominant thermal fre- 
quency at 300 K. The second-order terms qual- 
itatively show the same behavior as the first- 
order ones,— but there is a notable difference 
between the electric and the magnetic part: 
The second-order contribution to the magnetic 
LDOS at least is smaller than its first-order pro- 
cessor, although only by a factor which is not 
small compared to unity, whereas the magni- 
tude of the second-order electric contribution 
remains comparable to that of the first-order 
one, and even slightly exceeds it. At the bar's 
center, the first- and the second-order electric 
contributions amount to roughly 1.5 times the 
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Figure 2: (Color online) Zeroth-, first-, and second- 
order contributions to the electric (a) and to the 
magnetic (b) part of the LDOS above a gold sur- 
face structured with the nanobar (|39p as sketched 



in Fig. [TJ for 



at an observation height 



of 10 nm above the base substrate plane, 
different scales. 



Note the 



zero-order term. This behavior is not acciden- 
tal; it can be quantitatively understood with 
the help of an elementary consideration: At suf- 
ficiently short distances, that is, for z—h not too 
large compared to the profile width w, the local 
geometry equals that of a flat surface, shifted 
by h against the base plane. Therefore, the 
electric LDOS D E at such a point r = [0, 0, z]' 
is approximately given by the LDOS D^, per- 
taining to a perfectly flat surface^ through the 
relation 



(42) 



Now the distance dependence of Df s in the near 
field is given by 



Df s ([0A *]*) 



1 



(43) 



the strong decay of this electric LDOS with the 
third power of the distance clearly aids the local 
approximation (I42p . Thus, one has 



D E ([0,0,zf) 



1 



{z-hf ■ 

Expanding in powers of h, this yields 



(44) 



D E ([0, 0, zf) ~ 1 + 3 4 + ^ + °( h3 ) > ( 45 ) 

r/*~j *y*- 
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allowing one to estimate the ratios of the con- 
tributions appearing in different orders: 



(a) 



L>e.(°) ~ z ' 



£)E,(2) 
£)E,(0) 



6^ ■ (46) 



With h/z = 0.5, which is the value consid- 
ered in Fig. H one obtains £)E,( 1)^/^,(0) = 
_D E '( 2 ) /D E >(°) = 1.5, in quite good agreement 
with the exact numerical result. Generalizing 
this argument to any order n, one observes 



£>E,(n+l) 
L)E,(n) 



3 /i 



1 2 



(47) 



Thus, while the series may still converge for any 
z > h, convergence at short distances would be 
rather slow; for h/z = 0.5 the magnitude of the 
leading successive contributions, normalized to 
the zeroth-order one, is 1 : 3/2 : 3/2 : 5/4 : 
15/16 : 21/32 : . . . .In this example, terminat- 
ing the perturbation series at the second order 
means that one collects only about half of the 
exact value. 

In the case of the magnetic LDOS, one has 



Df s ([0,0,zf)^- z 



(48) 



so that here a larger area of the sample's surface 
contributes to the LDOS than in the electric 
case (|43p , implying that a local approximation 
analogous to Eq. (|4"2"j) cannot be expected to 
work as well as before. Ignoring this restriction 
and performing the analysis nonetheless, one 
ends up with 



£)H,(n+l) 
£)H,(n) 



h 

z 



(49) 



which, in view of the shaky foundation of the 
reasoning, still works satisfactorily, capturing 
both the correct trends and the orders of mag- 
nitude read off from Fig. [2] 

It is evident that these general findings do 
not depend on the specific dielectric properties 
of the material. Indeed, when considering a po- 
lar sample with a permittivity described by the 
Reststrahlcn formula 4 ^ 



(50) 



e{u)) = eoo 1 + 



and inserting parameters appropriate for gal- 
lium nitride 4 ^, namely, = 5.35 for the high- 
frequency permittivity, cl>l = 1-41 x 10 14 s _1 
and wt = 1-06 x 10 14 s _1 for the frequencies of 
the longitudinal and transversal phonons, and 
7 = 1.51 x 10 12 s _1 for the relaxation rate, and 
again taking the profile (|3"9"]l . we obtain Fig. 
which shows the same qualitative features as 
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Figure 3: (Color online) As Fig. [5] but for a sample 
consisting of gallium nitride (GaN). 



the previous Fig. [5] for the gold sample, al- 
though, of course, the scales are quite different; 
now the electric contribution dominates. Like- 
wise, the results do not seem to depend sensi- 
tively on the precise form of the structure: The 
corresponding data obtained for the reference 
profile (|40[) arc remarkably similar to the previ- 
ous ones, as shown in Fig. 21 and again confirm 
the simple estimates ([4"T|) and (|4T)|). For this 
reason, we only consider the gold nanobar (|39p 
in the following. 

To conclude the above discussion, in the 
cases studied so far the restriction to second- 
order perturbation theory already seems ques- 
tionable, although the strong dominance of the 
magnetic LDOS for metallic samples might still 
mask the problem with the electric one. This is 
potentially important for NSThM-applications, 
where typical probe-sample distances range 
down to a few nanometers. On the other hand, 
low-order perturbation theory may be expected 
to work reliably when the profile height h 
clearly is the smallest length scale of the prob- 
lem; according to the above reasoning, it should 
become better when increasing the observation 
distance z. In order to estimate the small- 
est z at which second-order perturbation theory 
might still give quantitatively good results for 
our model profile we plot in Fig. the 

ratios of the various contributions, evaluated 
above the bar's center at varying distance; here 
we also consider the regime z—h^>w where the 
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Figure 4: (Color online) As Fig. for a gold sample 
with the reference profile Q40|). 



short-distance estimates (|4T|) and ()49[) may no 
longer be taken for granted. As a rough guide- 
line, one may accept the truncation of the per- 
turbation series at the second order if the ratio 
of the second-order contribution to the zcroth- 
order one drops below 10%, say. For the domi- 
nant magnetic part this criterion is satisfied for 
z > 20 nm, while the electric part then requires 
z > 41 nm. (When reducing the acceptance 
limit to 5%, one gets z > 26 nm and z > 56 nm, 
respectively.) The potential problem of slow 
convergence here is expressed by the fact that 
the ratio of second- to first-order contribution 
decreases only rather slowly with increasing z. 
Still, in view of the relative smallness of the 
electric LDOS above metallic samples this does 
not seem to be essential. 



B. Resolution of an NSThM 

Scanning the LDOS at a constant height 
above the base plane, as done numerically in 
Figs. [3J [31 and HI corresponds to the constant 
height mode of operation of an NSThM. This 
is not an advantageous mode, for two practi- 
cal reasons: On the one hand, it is hard to 
realize with sufficient accuracy, on the other, 
it contains the risk of a probe-sample collision 
when scanning a surface with an unknown to- 
pography, almost always resulting in irrepara- 
ble damage to the delicate sensor A much 



Figure 5: (Color online) Ratios of the zeroth-, first-, 
and second-order contributions to the electric (a) 
and to the mag netic (b) LDOS, for w = 10 14 s" 1 at 
an observation point x = above the center of the 
bar ([39]). 

more favorable mode avoiding these complica- 
tions is the constant transfer mode, meaning 
that the sensor height is continuously regulated 
such that the detected heat current remains 
constant when moving the sensor over the sur- 
face; the information about the sample's near 
field then is embodied in the recorded sensor 
height. Note that this latter mode differs from 
the constant distance mode, which uses addi- 
tional information on the local distance of the 
sensor from the structured surface (obtained 
by electron tunneling spectroscopy) in order to 
keep that distance constant. That constant dis- 
tance mode was employed experimentally by 
Kittel et al.-^ numerical first-order results per- 
taining to this mode can be found in Ref. [27|. In 
contrast, a major benefit of the constant trans- 
fer mode lies in the fact that it exclusively re- 
quires heat-transfer information, so that there 
is no need to retain the sensor's ordinary STM- 
capability. Moreover, it allows one to assess the 
resolution limit of the NSThM in a simple man- 
ner. 

For modeling this constant transfer mode, 
we select some appropriate fixed value of the 
LDOS, and then calculate that observation dis- 
tance a at which this LDOS-value is reached. 
Only the sum of all contributions is of interest 
now. For consistency, we also require that the 
second- to zeroth-order ratio remains less than 
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Figure 6: (Color online) (a) Observation dis- 
tance above the one-bar gold sample described by 
Eq. (|39p (1 bar), and above the sample with two 
parallel bars specified by Eq. (|5ip (2 bars), deter- 
mined such that the second-order LDOS for ui = 
10 14 s _1 constantly keeps that value which is at- 
tained for the distance a g — 20 nm far away from 
the bars, (b) Definition of the quantities d and Ad 
employed for discussing the resolution of an ideal 
NSThM. 

10% for the magnetic LDOS, as the electric 
one does not contribute significantly here. Dis- 
cussing the NSThM's resolution power first re- 
quires the specification of a norm structure con- 
taining the length scale s to be resolved. Here 
we take two parallel gold bars of the form (|39f , 
as described by the profile 



J y ' |_exp [Ci(|x + 0.5s| - 0.5wi)] + 1 

(51) 

1 

+ e X p[C 2 (|x-0.5s|-0.5w 2 )] +1 ' 

the length scale in question is their separa- 
tion s. For our matter-of-principle calcula- 
tions we again choose h = 5 nm, together with 
W\ = u>2 = 30 nm and (i = C2 = 10 9 m . 

In Fig. [B^a) we display second-order results 
for both the one-bar geometry, and for the two- 
bar structure with bar separation s — 50 nm. 
Here the observation distance a is computed 
such that the LDOS remains fixed at the value 
attained for the distance a c ff = 20 nm above the 








separation distance s (nm) 

Figure 7: (Color online) Ratio Ad/d for three dif- 
ferent values of a e s (20 nm, 30 nm, and 40 nm), as 
functions of the bar separation s. 



base plane at positions far away from the bars, 
always assuming u> = 10 14 s _1 . 

In order to discuss the resolution of an ideal- 
ized NSThM, we make two further assumptions: 
First, we propose that the sensor is point-like, 
so that no effects due to the real sensor's ex- 
tension are considered, implying that we aim 
at the sensor- independent best possible resolu- 
tion limit. In reality, the finite sensor size will 
lead to a lower resolution. Second, we assume 
that the signal recorded by the device is pro- 
portional to the LDOS at the dominant ther- 
mal frequency, which actually appears to be 
quite a good approximation for metallic sam- 
ples.— We then take the ratio Ad/d, where d 
is the maximum difference a — a c ff encountered 
above each of the two identical bars, and Ac? 
denotes the difference between that maximum 
distance and the minimum distance adopted be- 
tween the bars, as illustrated in Fig. [(Jb). We 
now stipulate that the two bars can be resolved 
if Ad/d > r; this number r characterizes the 
sensitivity of the respective experimental set- 
up. 

In Fig. [7] we plot the ratio Ad/d for effective 
distances a e ff = 20 nm, 30 nm, and 40 nm, as 
functions of the bar separation s. If we take 
r = 0.2 for the sake of discussion, the resolv- 
able minimum distances are s m i n ~ 50 nm for 
a e s = 20 nm, s m i n ~ 65 nm for a e s = 30 nm, 
and s m in ~ 80 nm for a e ff = 40 nm. We em- 
phasize that these figures serve to illustrate 
the basic principle and should not be taken 
at face value. The resolution achievable with 
an actual NSThM device will also depend on 
the type of surface structure under investiga- 
tion; a further limit will be imposed by the fi- 
nite sensor volume. The key message, however, 
stands out clearly: When scanning an isother- 
mal, nanostructured surface with a Near-Field 
Scanning Thermal Microscope, one is able to 
resolve structures with linear extensions which 
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may fall orders of magnitude below the scale set 
by the dominant thermal wave length.— 



IV. CONCLUSIONS 

We have demonstrated that a numerical eval- 
uation of Grcffct's perturbation series for the 
scattered electromagnetic field at rough sur- 
faces^ is routinely feasible up to second order. 
This allows one to evaluate the local density 
of states above surfaces with arbitrary profiles, 
thus lifting the restriction to the limited class 
of profiles which can be dealt with analytically. 

The convergence properties of this series 
seem to warrant further analysis. While one 
may reasonably guess that low-order perturba- 
tion theory should be sufficient when the pro- 
file height h is by far the smallest length scale of 
the problem, the slow decrease of the successive 
contributions to the electric LDOS depicted in 
Fig. [5][a) for the smoothed, but still quite steep 
metallic model structure sketched in Fig. [JJ to- 
gether with the elementary estimates based on 
Eq. (|42p . are warning signs. While our results 
have been obtained for specific model profiles, 
it would be quite helpful to have mathemati- 
cally rigorous and sharp upper bounds on the 
higher-order contributions for any given surface 
structure. 

We did not discuss possible effects due to 
the nonlocal dielectric response of the sample, 



J 



which might come into play at distances of a 
few nanometers.— The question whether such 
effects would be detectable with a Near-Field 
Scanning Thermal Microscope (NSThM) de- 
serves further investigations. 

With respect to NSThM surface imaging, we 
have shown how to estimate the best possible 
resolution limit, attained for a point-like sensor. 
Here we have introduced a mode of operation 
characterized by constant heat transfer, giving 
access to isolines of the LDOS. It is now a major 
task to extend the preliminary studies reported 
in Ref. H2L considering surfaces with both regu- 
lar and random nanostructures, and to further 
explore the concept of near-field thermal imag- 
ing. 
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Appendix A: Calculational details 



In this Appendix we state the precise forms of the expressions which have been used in Subsec- 
tion [TlB] The vector introduced in Eq. ((3"4"]) , required for computing the first-order electric 
contribution (|3"3"|) . contains the products of the transmission coefficients i s (/c) = (2k z o)(k Z Q + fc z ) _1 
and t p (n) = (2^/ek z o)(ek z o + fc z ) _1 , together with a convenient prefactor: 



S (1) (k,k') = 



— X Q-i( k za+k' z0 )z 



kk' 167r 2 k z ok' 



zO 



t s (n)t s (K')K ■ K 1 

-* s («)*p(0^(« x K ') 



'pW'pfr') ^ - k z k' z {K ■ K 1 )] 



(Al) 



The other vector appearing in Eq. ([34| contains the traces of the dyadic products, 



A\ x \k, k ) — 



2„/2 



K ■ K 



(A2) 
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For computing the magnetic contribution, this vector has to be replaced by 



A tt (k, k ) — - 

A" A" ' 



KK 



A-(k 2 k' 2 - k z0 k' z0 (K, ■ k')) 



^■(k x k') 

In \ / 

^■(k x k') 



(A3) 



(2) 

The vector S\ {k, k') determining the second-order term (|35p is given by 



s¥\k,k') =i 



fc 2 (e - 1) e-^+^o)* 
167t 2 kk' k z ok' z0 



t s (K)t s (K')(k z + k' z )(K ■ O 
-t s (K)t p (K')^{k z + k' z )(K X k')* 

-* p («)4(k')^(^^ + n 2 fc 2 )(K x 



i p («;)ip(K')- 



hi 



/cV 2 (n% + k z ) 



(A4) 



whereas the vector (k, k', k") entering into the expression (|38p . and thus into the other second- 
order contribution ([3"T|). takes the form 



8tt 2 / 



i a (K)t B (K") 



■ K ^ly°y (« X «%(«' X K,'% 



-t s (K)tp(K") 



K'^+fc' fc: ^ x K > z 



-tp(K)i B (K") 



nfco n'' 2 -\-k' !, n k' 



tp(K)t p (#c") 



|tfr(K x X 



CA5) 
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